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ABSTRACT 

Recently, dense gases have been investigated for many engineering applications such as for turbo- 
machinery and wind tunnels. Supersonic nozzle design for these gases is complicated by their nonclassical 
behavior in the transonic flow regime. In this paper s method of characteristics (MOC) is developed for 
two-dimensional (planar) and, primarily, axisymmetric flow of a van der Waals gas. Using a straight sonic 
line assumption, a centered expansion is used to generate an invisdd wall contour of minimum length. The 
van der Waals results are compared to previous perfect gas results to show the real gas effects on the flow 
properties and inviscid wall contours. 


INTRODUCTION 

A m in i mum length nozzle (MLN) produces a uniform supersonic flow with a minimum ratio of 
throat height or radius to total nozzle length. Argrow and Emanuel (refs. 1, 2) present detailed discussion 
of MLN designs, comparisons of different types, computational flow field analysis, and engineering ap- 
plications. Until this paper, all MLN analyses, that we are aware of, have focused only on perfect gas re- 
sults. Several authors (refs. 3-10) have recently investigated nozzle flows for nonideal gases. References 
3-5 and 8-11 discuss steady flows of dense gases. Cramer (ref. 11) refers to these gases as BZT (Bethe- 
Zel’dovich- Thompson) to recognize the individuals that first theorized their nonclassical behavior. Such 
behavior includes local a minimum of the Mach number during steady isentropic expansion, expansion 
shocks, increase in the critical Mach number, and other nonclassical behaviors. BZT fluids typically have 
large polyatomic molecules with comparatively large specific heats. 

All the previously mentioned references that investigate BZT fluids discuss the thermodynamic 
condition that governs the classical or nonclassical behavior. The governing thermodynamic parameter is 



is die speed of sound, g, p, and I are the density, pressure, and entropy, respectively. The overbar indicates 
dimensional quantities. The parameter I' is referred to by Thompson (ref. 13) as the fundamental derivative 
of gas dynamics. 

For most fluids T > 0 under normal conditions, but for BZT fluids may have r < 0. The T < 0 
region of a BZT gas (dense gas) occurs in the dense gas region near the saturated vapor curve in the pv 
plane. Practical uses for the nonclassical behavior of BZT fluids include turbomachinery and heavy -gas 
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wind tunnels. Discussion of the various applications that may capitalize on the nonclassical behavior of 
BZT fluids can be found in refs. 3, 5, and 8-12. 

The most apparent application of a MLN designed for dense-gas flow is for heavy-gas wind tun- 
nels. Anderson (ref. 9) shows Navier-Stokes calculations for the flow of sulfur hexaflouride (HFg) over a 
NACA 0012 airfoil. This study investigated the feasibility of using this large-molecule gas for wind tunnel 
applications. While most of the previous references speculate that nozzles can be designed to produce su- 
personic flow of a BZT fluid, the calculations have been limited to one-dimensional cases (refs. 3-5, 8, 10, 
11 ). 


THERMODYNAMIC MODEL 
The van der Waals equation of state is 


RT a 
~ 7 2 


7-6 ~ 2 


where ?? is the specific gas constant, v is the specific volume, and 


8p c 64 P c 


The c subscript refers to conditions at the critical point. The following thermodynamic development and 
nondimensionalization scheme follows that of ref. 12. The enthalpy 5 and speed of sound a are given by 




ar - 
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Here e r is an arbitrary reference energy and 


<5 a /?/F v 


Because there is minimal temperature variation for the flows investigated, we assume a constant c v = c va 0 , 
where c va >is the ideal gas specific heat. The temperature variation for isentropic flow of a van der Waals 
gas is given by 



where the 0 subscript refers to stagnation conditions. 
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The nondimensional (reduced variable) form of the van der Waals equation is written as 


P - 


87 3 

3v - 1 v 2 ’ 


and the critical compressibility is 


Z c = 



with the reduced variables 





The nondimensional form of the other thermodynamics relations are 

3 


where 


r(1 + - % 


M = 


J2(h 0 - h) 


U - h - e r _ Z 

h ~WT' a 3 


m' 


b = A = J_ a * _a_ = 27_l 

8 Zc' a wr c V c 64 Z c 


METHOD OF CHARACTERISTICS FOR REAL GASES 

The method of characteristics (MOQ used by Argrow and Emanuel (ref. 1) assumes an isen tropic, 
irrotational flow of a perfect gas. For this case, the governing two-dimensional partial differential equation 
reduces to a set of four algebraic equations, two characteristic and two compatibility equations. For the 
axisymmetric case, the characteristic and compatibility equations form a set of four ordinary differential 
equations that are solved simultaneously. Details of the solution procedures can be found in ref. 1 . 

A MOC for the isentropic two-dimensional or axisymmetric flow of a real gas (ref. 14) is used for 
the present study. The method is completely general and we use the van der Waals thermodynamic model 
more for simplicity than for accuracy. The governing partial differential equations (PDEs) are given by the 
gas dynamic equation, the irrotationality condition, and the speed of sound relation, 

{Vl - a 1 )V x + (V 2 - a 2 )Vy + 2V x V y ^ - = 0, 
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dVy 

dy dx ’ 
a — a(V) = a(V,, V y ), 

where a = 0 for two-dimensional flow or a = 0 for axisymmetric flow, x and y are the axial and radial 
(transverse) coordinates nondimensionalized with respect to throat radius or half-height, V x and are the 
corresponding velocity components, and V is the velocity magnitude, all velocities are nondimensionalized 
in the same manner as the speed of sound a. Along the characteristic lines, that correspond to Mach lines in 
the flow, this system of PDEs reduces to two ordinary differential equations called compatibility equations. 

Fig- 1 is a schematic of the supersonic flow field showing the MLN geometry and the flow geometry 
associated with an arbitrary point on a streamline. The angle 9* is the initial inclination of the supersonic 
contour. The left and right running characteristics designated as C + and C_ are Mach lines inclined at the 
Mach angle fi with respect to the velocity vector V. The corresponding characteristic and compatibility 
equations are 

C ± characteristic equation: = A ± = tan(0 ± /*), (la) 

C ± compatibility equation : (V| - a 2 )dV x± + [2V x V y - (V 2 - a 2 )X ± ]dV y± - a^-dx ± = 0 . (16) 

The gas is assumed to enter the supersonic portion of the nozzle along the straight sonic line OA 
uniform and parallel to the axis. It is then expanded and accelerated (or possibly decelerated in the r < 0 
region) through the nozzle and exits the nozzle with a uniform flow crossing the terminating characteristic 
BC at the exit Mach number Mf . For the two-dimensional nozzle, the centered expansion generated by the 
sharp throat is a Prandtl-Meyer expansion. In the axisymmetric case, the flow at the wall is locally two-di- 
mensional, thus at the throat the expansion is locally Prandtl-Meyer. Construction of the flow field begins 
by discretizing the centered expansion into equally-spaced velocity increments. Each velocity increment 
A V has an associated isentropic flow turn angle increment AO: The angle increment A 6 is computed from 
the relation 


0 2 - d x 


i: 


l f\p- - 1 


dV . 


( 2 ) 


For a perfect gas, the AM associated with an angular deflection A 9 can be easily determined from the 
Prandtl-Meyer function. For a real gas, the Prandtl-Meyer computation requires the solution of a system of 
ordinary differential equations as shown by Cramer (ref. 10). We avoid the Prandtl-Meyer computation by 
using the angle-velocity relation, equation (2). 

With the position of the throat specified and the velocity components (V x ,Vy) computed from V and 
9, the necessary independent variables are determined at the throat. Equations (1) are solved simultaneous- 
ly using the second-order average-property Euler predictor-corrector scheme described in ref. 14. The 
characteristic net of the kernel region OAB and the transition region ABC is constructed using the unit pro- 
cesses described in ref. 1. The wall contour corresponds to the streamline tha t passes through points A and 
C. This streamline is determined by also using the average-property Euler predictor-corrector scheme to 
integrate the equation 

If = tan9 '" 

from the initial condition, 9 W = 9* at * = 0 to the exit condition 9 W = 0at.y = 0. 
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The accuracy of the characteristic computations and the wall contour are affected by the spacing of 
the characteristic nodes. The C_ characteristics that emanate from the comer at point A reflect from the axis 
as C+ characteristics and bend away from the center of the expansion. This causes relatively large spacing 
between characteristic node points and the subsequent computed wall points just downstream of the throat 
where the gradients are largest. To alleviate this situation, the characteristic net is compressed toward the 
sonic line by decreasing the step size of the speed increments in the centered expansion discretization. The 
Prandtl-Meyer and turn angles in ref. 1 are discretized and compressed in a similar manner. 

RESULTS 

For the nominal case of To =1.01, <5 = 0.02, and vq = 0.70, dVadV - 10 -6 for equation (2). The 
first characteristic is chosen at dV, then the step size is increased until a characteristic is generated at every 
206 V. This allows characteristics to be compressed towards the sonic line for improved accuracy. For the 
nominal case, the limits of integration for equation (2) are set at the sonic speed V\ = 0.26 to V 2 - 0.9^,^ 
with Vmax = 9.93. Characteristic spacing in the transition region is controlled by an aspect ratio that keeps 
the shape of the characteristic cells as uniform as needed. For the relatively short nozzles produced in this 
study, the aspect ratio was set at 0.5. Reference 1 gives a complete description of how the characteristic 
compression and transition region aspect ratio affect the overall computational accuracy. 

Figures 2 show the wall contour and the variation of M, q, and T along the axis and along the wall 
for a two-dimensional nozzle with To =1.01, 6 =0.02, and vq =0.70, 0.85. Figure 2(a) indicates the increase 
in nozzle length required for the more dense gas. Note that 6* is fixed at 2.5° for the comparisons in Figs. 
2. This is near the maximum value of 0* where the C. characteristic, near AB of Fig. 1, begin to overlap 
other C. characteristics in the kernel, or there is overlapping from C+ characteristics near BC with both 
instances producing an oblique expansion shock. Once a shock occurs (compression or expansion), the 
flow is no longer isentropic and the MOC cannot be applied without some special procedure to account for 
the placement of the shock (We do not to incorporate such a procedure in the MOC used for this study.) 

For Figs. 2(b-d), the axis curves are terminated at the end of the kernel because the uniform flow 
region begins at point B as shown in Fig. 1. The wall curves extend to the end of the contour at point C as 
also shown in Fig. 1. The density q decreases smoothly along the axis and wall for both cases as shown in 
Kg. 2(b). Figure 2c shows that the gas is expanded into a /'< 0 region that extends through the nozzle exit 
for the nominal vq =0.70 case. This produces a uniform supersonic flow with T < 0. For this case, the 
Mach number reaches a maximum of about 1 .88 along the axis and along the wall, as shown in fig. 2d. 
This agrees with the quasi one-dimensional results obtained by Cramer and Best (ref. 12). 

A comparison of the wall contours at the maximum Mf for the two-dimensional and axisymmetric 
nominal cases is shown in Fig. 3. Note that in the contour plots, the y-axis is not to scale. Although diffi- 
cult to see, the axisymmetric contour contains an inflection point just downstream of the throat that is not 
present in the two-dimensional contour. Also, Q* is larger for the two-dimensional case as will be shown in 
more detail in the following figures. The results of this plot agree with the perfect gas results of Argrow and 
Emanuel (ref. 1). 

Figures 4(a) and 4(b) show the 0 vs. Mf variation for the axisymmetric and two-dimensional cases, 
respectively. These figures show the nominal case, a perfect gas with nominal conditions, and nominal 
conditions with vq = 0.85. Varying stagnation conditions for the perfect gas does not effect the nondimen- 
sional results. The dense gas case, vq = 0.7, shows BZT gas behavior, reaching a maximum Mf of about 
1.88 for^ 6 near 1° for the axisymmetric case shown in Fig 4(a). Then Mf decreases until reaching a maxi- 
mum 6 value of about 2.5° . The maximum Mf for the two-dimensional case is at 6 * a 2.5° . The decrease 
of Mf for 6 > 2.5° is not shown in Fig. 4(b), although the maximum 0*for this case is about 3.5°. 

Figures 5 and 6 also show 6 vs. Mf for the for the nominal axisymmetric case. Figure 5 shows the 
effect of varying 7*0, keeping other nominal conditions fixed, compared to a perfect gas. This shows that 
slightly increasing Tq moves the flow away from the dense gas region. Figure 6 shows the effect of varying 
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c v /R, keeping the other nominal conditions fixed, compared to a perfect gas. This shows a tendency to ap- 
proach petfect gas behavior as cJR is decreased. 

Figures 7 and 8 show the variation of the kernel length xg and nozzle length xf vs. Mf, respectively. 
Note that for the two-dimensional case, the figures show the increase to the maximum lengths but the de- 
crease as 6 increases (and Mf decreases from the maximum) is not shown in these plots. The BZT gas 
behavior is indicated in the nominal cases. Because the gas is more dense, the lengths are longer and reach 
a maximum before decreasing. 


CONCLUSIONS 

The MOC is applied to the steady isentropic flow of a dense gas. A method is presented for gener- 
ating in viscid MLN contours to produce a uniform supersonic flow. The MLN procedure presented is Um- 
i' id to cases where there is only one sonic point Inclusion of more than one sonic point will require the 
coupling of the present method with a subsonic contour design procedure. In order to continuously isen- 
tropically expand the gas from a stagnation state to M —*■<*> t reference 12 shows that for the nom inal case 
investigated in the present study, M first reaches a local maximum value of about 1.88 before decreasing to 
a supersonic minimum of about 1.05. We have shown that it is not possible for a nozzle to use a single 
centered expansion to accomplish this. At best, the single centered-expansion MLN can generate an isen- 
tropic expansion slightly beyond the local maximum before an expansion shock is generated. The maxi- 
mum allowable 6 appears to occur as the r < 0 region approaches the throat, resulting in a shock. We 
speculate that a design that employs two centered expansions separated by a finite converging wall may be 
able to accomplish the expansion in a minimum length. This will be investigated in the future. 

It was shown that the nozzle designs presented may be used to produce a steady supersonic flow of 
a BZT gas in a r < 0 state. This raises the possibility of producing interesting wind tunnel experiments to 
study the supersonic external flow of a dense gas over airfoils and other aerodynamic shapes. 
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B 


Fig. 1 MLN supersonic flow field geometry, upper half-plane. 
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Fig. 2 (a) Wall contour, (b) density, q> 
tions; (d) Mach number. A/, distribution! 


case, nominal and for vq = 0 . 85 . 


V, = 0.70, 

v, - 0.85 


(a) 
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istributions (c) fundamental derivative, r, distribu- 
along the wall and the axis, all for the axisymmetric 















Fig. 3 Comparison of two-dimensional and axisymmetric contours for maximum 
Mf , nominal case. 
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^8* ^ Kernel length xg vs. Mf, nominal two-dimensional and axisymmetric 
cases with variable vq. 



Fig. 8 Nozzle lengthy vs. M f , nominal two-dimensional and axisymmet- 
nc cases with variable vq. 


341 




itinwimriffn 




